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SUMMARY 


IXiring  the  first  six  months  of  the  subject  contract,  the  puiposes  of 
the  research  protram  wore  to  Investigate  properties,  of  an  earthquake  source, 
to  provide  a  phy:  leal  basis  for  discrimination  between  different  typer,  of 
seismic  events.  The  two  main  approaches  used  were: 

1.  To  compute  theoretical  seismograms  which  fit  observed  ground  motion 
in  the  near-source  region  of  an  earthquake. 

To  investigate  shear  dislocations  in  the  context  of  crack  theory, 
requiring  tliat  a  sociated  shear  stresses  on  the  moving  fault  surface  be 
dynamically  satisfactory  in  terms  of  the  mechanics  of  shear  failure . 

Research  on  the  first  of  these  approaches  is  heavily  concerned  with 
computational  methods.  The  required  theoretical  seismograms  involve,  at 
each  station,  and  at  ea.ch  time,  a  triple  integration.  Veiy  considerable 
savings  In  computational  effort  can  be  made,  if  maximal  efficiency  Is  attained. 
By  carefully  arranging  the  sequence  of  Integrations,  and  by  using  a  valuable 
number  of  points  for  each  integration,  we  have  achieved  a  tenfold  improvement 
in  the  speed  of  these  computations.  We  have  also  been  able  to  extensively 
generalize  the  types  of  source  motions  for  which  resulting  near-field  ground 
motions  may  be  computed.  Present  computations  of  this  type  are  still  limited 
by  the  assumption  that  the  earthquake  is  talcing  place  within  an  infinite  elastic 
medium,  rather  than  near  the  free  surface  of  an  elastic  half  space.  Our*  present 
efforts  are  directed  towards  dropping  this  assumption,  thus  allowing  the  effect 
of  the  free  surface  to  bo  computed.  We  have  shown  that  this  will  be  important, 
because  we  have  found  that  elastic  surface  waves  do  contribute  significantly  to 
the  strong  pawl-motion  in  the  near  field  of  an  earthquake. 


Out'  of’  tlto  least  well  known  prirvui.  tops  ol’  an  o:u *11  iquake  }:;  11,.-,  stiver. 


drop.  Theoretical  cons id  tions  li uli « vit-c'  that  1.1 1 1  :•  quantity  Is  proportional 
to  the  pnrticle  ve Loot ty  at  the  souiro ,  and  the  coir.lant  of  pi*oportlonulity 
Inc  been  given  In  the  Literatim'  ac  a  complicated  formula  involving  the  r up- 
tun'  velocity  for  the  earthquake,  an  well  as  the  ela:  ic  parameters  of  the 
surrounding  medium.  We  hav  ’  :  uccoed*  d  in  evaluating  thin  relationship  for 
a  variety  cf  nipt ur-e  velocities. 

fair  piii din-;  principle  behind  our  second  major  research  approach,  the 
theoretical  study  of  shear  dislocation.",  which  are  conipntnble  with  the  mechanics 
of  shear  failure,  is  as  follows.  The  hr  vledge  of  the  physics  of  failure  on  a 
fault  surface  gives.  us  insist  it  into  the  static  ef  stress,  and  its  time  hi  "tory 
during  niptuie,  in  the  source  region.  However,  to  calculate  local  ground 
notions,  we  must  first  infer  the  dislocations  on  the  fault  surface  which 
arise  from  the  state  of  stress.  This  deduction,  of  the  dislocations  frem  tlie 
stress,  is.  a  mathematical  problem  .in  tlie  fielu  of  singula:-  integral  equations. 

It  is  a  subject  about  which,  for  three  dimensional  problems,  little  is  lenown. 
Although  the  actual  elasticity  problem  we  here  v  ok  to  solve  is  one  which 
involves  vectors,  we  ha r  fortunately  found  an  analog  scalar  problem  which 
well  illustrates  the  essence  of  difficulties  in  the  theory.  Our  main  research 
effort  heiv  has  been  directed  towai-u.  finding  what  sub-region  of  the  fault 
surl’ace  is  controlling ,  via  its  stress  history,  the  dislocation  at  each  point 
on  the  fault.  Although  our  in:. in  results  to  date  aiv  analytic,  there’  is  every 
indication  that,  the  problem  wi  11  he  aiiien-.ble  to  a  joint  computational-analytic 
approach.  Certainly,  it  appears,  that,  computers  of  tlie  pwsent.  general  Ion  are 
nowhere  near  big  enough,  for  this  problem  to  be  solved  by  numerical  studies  alone. 


3. 


K, .port; t •  on  CcmpulhUoiinl  rrocoduros.,  for  (V nonl i Thonrot;^ ; ' «J 

from  a  {'Ivon  Fault  Surface  on  whirl  hi  slo'-d.ir,f.  J  «_A1 ready  Sped  fled  ar  a 

Function  of  fos i t ion  and  Time 

In  thin  section ,  wo  prevent  the  Improvements  in  recent  programing.  Thin 
work  stems  from  the  clastodyrnmic  ropres  tat  ion  theorem  in  the  foim  Civ''n 
by  dc  Ibop  (193-3).  A  form  of  this  theorem  appropriate  for  the  representation 
of  a  faulting  source  in  an  infinite  homogeneous  medium  in  di scussed  in  detail 
by  lhrV  ’ll  (i960) .  The  source  in  described  in  terms  of  a  shear  dislocation 
propagating  over  the  fault  surface.  Thin  dislocation  can  be  either  an  edge 
dislocation,  or  a  screw  dislocation,  and  either  of  these  two  types  will  of 
course  generate  a  thr  -component  vector  displacer  "nl  at  every  point  wit  u.n 
the  medium.  Tnun,  six  possible  scalar  equations  must  be  considered,  to  include 
all  poss  .uie  displacements  from  either  of  the  possible  dislocations.  An  example 
of  these  formulae  would  be  the  following: 


/V  6*  fr/B  ,  , 

u«(*,t)  =  (B2/M  Y  '6(5y  2  "  VW*  "  t  )t;  dt 

l  ~  JJS  3  ,  1  Jr/a112 


+  2(6y  2  -  1)  (ar)“2D  (f,  ,f,  ,t  -  r/a) 

1  112 

-  3(,Iy  2  -  1)  (3r)"2D  U  ,K  ,t  -  r/B) 

1  1  1  2 

+  ?Y  2  (a^r)-1D  (f.  ,K  ,t  -  r/a) 

1  112 

-  (?y ^ 2  -  1)  (B3J*)"1Dj(C1,C2,t  -  r/B^d^df^. 


In  this  equation, 

S  =  Fault  plane  area. 

u  -  (u  u  ,  u  )  =  Cartesian  components  of  displacement  measured  from  the 

♦V  1  ’  2*  3 

initial  state. 


'I. 


)'  -•  (x  ,  x  .  x  )  -  C'ii‘1' .  i'in  coordinate;  or  point.  at  which  u  i:  to  he  evaluate'! 

^  1  2 '  3  ~ 


C  ~(5,t  >  0  )  =  Cart'  dun  coordinates  of  point  of  integration  on  8. 
'X'  12  3 


p  =  density. 

ci  =  P-vnve  velocity. 

£  =  S-weve  velocity. 

n  =  (n  ,  n  ,  n  )  ~  unit  nor  rial  on  8. 

A'  1  2  3 

D  =  (D  ,  D  ,  D  )  =  di  sp]  uci.riont  ill  .".continuity  across  8. 

^  x  1  2  3 

'Hie  super*. ’-or  Inc fl ,  on  tlie  displacer'  >nt  calculate!  by  equation  (1),  indicates, 
that  the  assumed  dislocation  on  S  is  an  edge  dislocation,  for  which  JD  =  (D  ,  0, 

0)  only. 

Equation  (1),  together  with  the  fiv  other  similar  equations  for  the  other 
displacements  component j  ^ three  with  the  screw  dislocation),  has  been  very 
widely  used  for  computations.  It  has  almost  always  been  assumed  that  D  has  the 

A/ 

form  of  a  rrurp  function,  the  dislocation  at  any  point  of  the  fault  surface  thus 
growing  steadily  with  time  (after  the  rupture  front  has  arrived),  until  the  final 
fault  offset  is  reached,  at  which  time  the  fault  is  assumed  to  loch. .  With  this 
special  assumption  for  tire  time  history  of  the  dislocation,  the  basic  computing 
formulae  nay  Ire  further  reduced  from  the  form  (1),  becoming  even  more  amenable 
to  programming;.  However,  equation  (1)  does,  contain  the  generality  permitirig  use 
of  possibly  more  realistic  dislocations .  It  will  be  noted  that  two  of  the  terms  4. 

in  this  basic  equation  contains  the  time  derivative  of  the  dislocation.  Differentia¬ 
tion  is  an  operation  which,  when  done  numerically ,  is  liable  to  be  performed  in¬ 
accurately.  Indeed ,  for*  some  of  the  dislocations  jusldcral  below,  the  time  deri¬ 
vative  actually  contains  singularities  (which,  however,  ;u*e  integrable) .  The 
numerical  work  requires  that  the  resulting  derivative  be  integrated  over  the 


taull  area,  and  thin  nan  ho  regard<  cl  a.’  a  smoothing  process.  However,  we 
kava  greatly  .Improved  the  slab •  J j ty  ol‘  oar  cunpulat.icns  hy  reversing  the  order 
o'  th(  (L  fferentiation  an  I  Intoj  nation.  A  i '(.'levant  term  In  equation  (1)  jo 


time,  treated  an 


constant  time; 


•  <1  i  2y  2  (u 
at  Jj  i 


3i’)”1D  (f,  ,  t  ,t  -  r/u)dr,  d£ 
112  1  ? 


Hie ^reversal  of  operation;;,  as  net  out  in  equation  ?,  has  yielded  fourfold 
improvement  over  conventional  prograsrriing  of  equation  1.  That  is,  only  one 
fourth  the  number  of  points  on  1  are  needed,  to  give  answer’s  with  similar 
accuracy.  Hero.  importantly,  this  order  of  nan *ri  'al  procedur -s  (integration 
first,  then  differentiation)  docs  permit  the  use  of  dislocations  with  more 
realistic  time  histories,  and  for  which  there  may  numerically  be  a  singularity 
in  the  particle  velocity  of  the  fault  surface. 


In  the  convent Ioaal  progr, aiming  of  equation  1  (see,  for  example,  Anderson; 

197^0,  a  fixed  number  of  points  ;  e  taken  in  the  £  and  £  directions  on  the  fault 

1  2 

surfaces.  For  long,  thin  faults,  perhaps  as  few  as  four  points  are  taken  in  the 
direction.  Since  the  basic  computations  (involving  a  total  of  3  integrals: 
see  equation  1)  have  to  be  performed  for  each  moment  in  time,  it  would  clearly 
be  advantageous  to  allow  a  different  spacing  to  bo  used  for  the  points  on  tlie 
fault,  for-  the  calculation  at  each  moment  of  time  at  the  receiving  position. 

To  accomplish  this,  we  now  use  Romberg  integration  for*  each  of  the  tliree  Integrals. 
(Romberg  integration  uses  successive  doubling  of  the  number  of  intervals  used  in  a 
Fr' veil  Integration  problem.  Hie  result  for  each  choice  of  intervals  is  stored.  As 
the  number  of  intervals  is  increased,  a  sequence  of  approximations  to  tin-  desired 


(». 


Integral; ion  ;i a  found.  Of  course,  the  do  ind  result  would  be  obtained  in 
the  Unit  an  intoi-val  sine  lend;  to  0.  From  the  sequence  of  approximations 
to  a  drained  Integral,  the  algorithm  predict:.;  the'  limit  to  which  the  approxi¬ 
mations  are  tending.) 

The  improvements  I'eeulting  from  use  of  triple  Romberg  Integration,  together 
with  the  sequence  of  operation  ;  stated  In  equation  ?,  have  resulted  in  a  tenfold 
improvement  in  computation  time  for  a  given  problem  (at  a  stated  level  of  accuracy). 
Such  a  saving  is.  highly  significant,  since,  once  a  program  of  this  tyne  has 
successfully  been  debugged,  it.  may  well  be  used  for  problems  requiring  major 
amounts  of  computing  time. 


/ 

* 

A  sevei'e  check  on  the  accuracy  cf  our  present  program  i;  afforded  by  compari¬ 
son  with  a  particular  problem  already  described  by  Richards  (1973a,  b).  These 
papers  describe  a  special  method  for  finding  theoretical  seismograms  which  result 
from  a  growing  plane  elliptical  shear  crack. Oar  new,  and  much  more  general,  method 
for  doing  these  computations  does  successfully  duplicate  tip  wave  s tapes  previously 
obtained  by  the  specialised  procedure. 


The  main  weakro  so  of  our  present  computational  method  is  its  underlying 
assumption  that  the  f;  siting  occurs  in  an  infinite  h  logeneous  elastic  medium. 
This  rules  t  any  effects  due  to  the  Earth's  free  surface.  For  example,  surface 
waves  are  excluded,  and  P-5V  coupling  (in  body  waves)is  also  excluded.  In  fact, 
the  energy  radiated  from  a  rupturing  fault  may  well  couple  efficiently  into 
surface  waves.  Wo  justify  this  speculation  as  follows: 


Con:1. idor  tin*  nit. u;;f  n  dcpiri.nl  in  Figure  L.  Thu-,  shows  pru  uiietorr.  i'or 
a  propagating  souiv. ■  within  an  elastic  medium,  an!  general] ng  surface  waver, 
'ilie  displacement  rerponre  ai.  the  receiver  may  be  i -egarded  r .imply  an  a  super¬ 
position  of  all  frequency  component r. ,  and  all  positions  occupied  by  tiie  propa- 
patlnp  r.ource.  Thus, 


A 


u(x,  l)  = 


J 


A  (t ,  x,  corO  ,  03)0^'’  d^dw 


(3) 


in  which* the  phase  factor  l  in 
*  =  w!  t  -  1  dn  -  kX(f„  x). 

J0  vTny 


-{ 


(4) 


Note  that  ‘there  equations  permit  a  rupture  velocity  V  which  may  vary 
along  the  path  of  faulting.  In  evaluating  the  double  integral  of  equation 
3,  it  may  be  expoc  ed  that  the  major  contribution  will  come  from  the  point 
of  stationary  phase,  obtained  for  each  time  t  by  equating  partial  derivatives 
of  4>  (with  respect  to  C  and  u)  to  0.  Thus ,  the  major  contribution  is  from 
values  of  to  and  S  such  that 


t  - 


'o 


dn 

YTn) 


X 

u 


and  u  =  V(U  F  (S,  x,  cosO) 
k 

where  F  a  -  DX  -  /  1  +  x?sin20 
DC  ^  (xcorO-02 


(5) 

(6) 
(7) 


mid  ’J  is  group  velocity. 


Tito  function  F  has  valu'  approximate ly  < -cju* 1 1  to  1  Jl1  tl to  ancle  0 


s  small.  For  distant  receiver.  ,  F  tent  to  the  valu  •  coot. 


Equation  5  states  that  surface  wave:;  ;uvive  at  a  time  given  by  the 


time  to  rupture  to  position  £,  plus  the  tine  taken  to  travel  distance  X 


with  group  velocity  U.  Equation  6  states  that  the  frequency  of  such  surface 


waves  is  tint  for  which  the  phase  velocity  equals  rupture  velocity  times  the 


function  F,  which  has  a  value  between  0  and  1.  These  results  can  be  summarised 


by  saying  that  coupling  into  surface  waves  will  be  efficient  if  rupture  velocity 


is  somewhat  greater  than  the  surface  wave  phase  velocity.  Present  estimates 


of  tin  rupture  velocity  of  faulting,  ah  t  2  or  3  kilometers  per  second,  are 


indeed  about  the  same  as  values  for  surface  wave  phase  velocities. 


■ 


Ccr»nistf>nl  v.’i  Ui  a  given  Time  History  of  Ht  rear 


rn ie  most  commonly  used  dislocation  in  earthquake  source  theory  probably 
is  tint  due  to  Haskell  (1969) ,tbe  "ra  :p  function"  described  in  tin  section 
above.  However,  tills  time  history  probably  would  result  in  unacceptable 
stress  singularities  on  the  moving  fault  surface  itself.  A  dislocation 
function, is  lsicwn,  for  which  tlie  associated  stresses  are  constant  on  the 
moving  part  of  the  fault,  arid  this  function  Ins  been  described  by  Kostrov 
(1G6;1),  and  Richards  (1973a,  b).  However,  tills  dislocation  is  relevant  only 
to  die  beginning  rupture  notions,  and  docs  not  describe  the  way  in  which  fault 
motions  cense.  As  described  in  the  original,  proposal  for  the  present  contract, 

S’ 

one  of  the  least  understood  aspects  of  earthquake  source  theory  is  the  way 
in  which  fault  motions  cease.  A  fully  numeric/1  approach  to  this  research 
problem  is  unlikely  to  be  successful,  since  three  spatial  dimensions  and  one 
tine  dimension  must  be  studied,  and  this  total  of  four  dime  ions  imposes 
enormous  memory  requirements  for  retaining  the  grid  on  whi  ii  finite  elements 
or  finite  difference  computations  are  based.  Frasier,  in  a  paper  presented 
at  the  197,1  Annual  Heeting  of  the  American  Geophysical  Union,  indicates  tint 
even  the  computer  ILL TAG  9  is  inadequate.  Thus ,  combined  numerical  and  analytic 
methods  provide  the  best  chance  for  a  successful  approach. 

Our  analytic  research,  into  the  relationship  between  stress  in  an  earth¬ 
quake  source  region,  and  the  resulting  d'slocation  across  the  fault  surface, 
has  been  advanced  by  finding  an  analog  problc  i  of  relative  simplicity.  Tills 
analog  problem  Involves  a  scalar  field, rattier  tlnn  the  vector  displacement  field 
of  elasticity.  The  role  of  traction  is  played  by  a  directional  derivative  of 
the  scalar.  We  turn  next  to  a  basic  exposition  of  Oils  analog  problem: 


.10. 


Suppose  Mat  4  -  <|‘(x,  t)  .in  seine  scalar  function  oL‘  space  and  time  winch 
satisfies  the  basic  wave  equation 
V2*  =  1  3^, 

v2  at7  (8) 

within  a  region  V,  and  yG(x)  is  the  initial  (t^O)  static  value  of  <p.  At 

time  t  =  0,  <p  begins  to  experience  a  growing  dislocation  across  an  internal 

surface  £ ,  which  nay  be  regarded  as  the  fault  surface,  and  is  contained  in  V. 

However,  -although  <p  is  discontinuous  across  £,3^  is  kept  continuous  (in 

an 

anaology  with  the  requirement  of  continuity  of  traction  across  a  fault 

surface).  The  problem  Is,  first,  to  find  the  resulting  4  (x,  t)  throughout 

V  in  terms  of  the  discontinuity,  in  <  on  I.  The  later  problem  will  be,  to 

specify  the  value  or  34  on  £,  and  to  deduce  an  associated  discontinuity  in 

3n 

4.  (This  second  problem  is  the  analogue  of  specifying  the  history  of  trac  ion 
on  a  rupturing  fault  surface,  and  finding  the  related  discontinuity  in  dis¬ 
placement  . ) 

An  important  role  in  these  problems  is  of  course  played  by  Green’s  function 
G(x,  t;  £,  x),  which  in  V  s:  isfies 

V2G  =  1  32G  +  S  (x-06(t-r),  (9) 

v2  at7 


the  S's  being  Dirac  delta  functions.  The  "infinite  space"  Green’s  function 
satisfying  (9)  is 


G  =  -  1  6  (t-x-R) 

PtR  v 


where  R  3  s  the  distance  |x-c|. 


V 


We  luvi  found  the  solution  lo  tire  I'.i  rt  problem  in  three  ditTeivnt;  fom;, 
eaeli  of  which  has  merit..  Tix.-y  ;u-o 


<!'  (x,t)  « 


Jj  cogs  f  i)  5(e,  ,hJ1  *1  U.  t.JO  / 

Y;  'inR  l  v  L  v  »■*  i;  L  ~  vJJ 

(>'  (■ 

-  1  jiioso  1  I  $(c,  t-R)  dz 
DK  L  R  i  r  vJ 


i  d): 
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in  which  [  ]  denote  the  discontinuity  taken  acn  Z  in  the  direction  of  n 

/.  * 

tire  latter  being  the  unit  normal,  and  0,  R,  x,  £  arc  as  shown  in  Pietro  2. 

Virtually  all  the  practical  problems  contained  in  the  elasticity  problem 
are  present  also  in  the  uses  of  equations  (11),  but  these  analogs  forms  do 
not  have  the  algebraic  complexity  of  the  vector  problems.  A  list  of  these 
problems,  now  being  actively  pursued,  include  the  following; 

(i)  io  r„udy  the  singularities  in  (11),  as  tire  receiver  position  x  is  taken 
ever  closer  to  the  fault  surface  E.  Tire  problem  here  is  tint  the  quantity  R 
tends  to  0  lor  some  part  of  the  surface  integration.  Since  the  integrand 
contains  tcras  like  R~2,  this  implies  non-integrable  singu lorities.  In  this 
connection,  it  is  important  to  note  that  if  the  term  in  square  brackets,  in  the 
last  ol  equations  11,  is  a  constant,  then  tire  resulting  integral  is  nothing  but 
the  solid  angle  subtended  by  Z  at  x.  'ihi s ,  of  course,  is  simply  Pi r,  if  x  is 

9  /'v 

on  tire  fault  surface . 

(ii)  To  study  the  region  or  fault  surface  actually  contributing  to  the  integrals 
at  a  specified  time.  The  point -hrre  is  that  all  the  integrands  in  equations  11 
(and  in  Haskell's  equation  (1)  above), arc  evaluated  for  retarded  time.  This 


retarded  time,  if  H  is  large  enough,  can  Mi-  n  be  no  early,  Unit  for  the 
point  t  or  integration,  1  h  >  rupture  front,  has  not  ye  t  arrived,  and  the 
dislocation  is  nero .  Thus,  a  large  part  of  the  fault  surface  may  not 
actually  he  contributing  to  the  Integration.  The  parte  of  Y.  which  contribute 
to  tiie  integration  are  defined  by  £  valu  a  satisfying  tl«  Inequality 

tsS 

t  _  R  .>•  t(£)  (12) 

—  **  fj 

V 

where  t  is  the  time  at  which  the  rupture  arrives  at  £ .  We  have  examined  the 

case  of  earthquake  motions  initiated  on  the  fault  surface  5:  by  a  dislocation 

which  begins  at  a  point,  and  grows  behind  a  circular  rupture  front,  which 

expands  v/ith  rupture  velocity  c<v.  The  sub-region  of  l  which  influences 

tiaaction  at  a  particular  point  on  l  can  then  be  feun.i  arnly ti cally ,  by  solving 

ecuation  (1?)  for  (£  ,£  ). 

1  2 
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Figure  1:  Parameters  for  a  propagating  point  source,  generating  surface  vnve-o 
at  the  receiver  position. 


